#include <iostream>
#include <opencv2/core/core.hpp>
#include <opencv2/highgui/highgui.hpp>
#include <opencv2/gpu/gpu.hpp>

#ifndef DDMCMC_Vector
#define DDMCMC_Vector

template<class Tp>
class Vector{
private:
	cv::Mat mData;
	int nSize;
	bool bRowVector;

public:
	// construction & destruction
	Vector(){};
	Vector(Vector<Tp>& a); //untested
	Vector(const int size){ SetRowVector(size); };
	~Vector(){};


	// operators
	Vector<Tp>& operator=(Vector<Tp>& x);
	Tp&			operator()(const long index);
	Tp&			operator()(const int row, const int col);

	Vector operator+(Vector<Tp> &x);
	Vector operator-(Vector<Tp> &x);
	//Vector operator*(const Matrix<Tp> &x);
	//Vector operator*(const Vector<Tp> &x);

	Vector		operator+(const Tp y);
	Vector		operator-(const Tp y);
	Vector		operator*(const Tp y);
	Vector		operator/(const Tp y);

	friend std::ostream &operator<< (std::ostream &strm, Vector& vec){ return strm << vec.getData() << endl; };
	//void operator<<(std::ofstream &strm);


	// operations
	inline long size(){ return nSize; };
	Tp Data(const long index);

	int		rows();
	int		cols();
	void	SetDimension(const int row, const int col);
	void	Zero();

	cv::Mat		getData();
	void	InitValue(Tp value);
	void	setData(cv::Mat mNewData);
	void	SetRowVector(const int size);
	void	SetColVector(const int size);
	bool	RowVector(void);
	bool	ColVector(void);
	bool	ExpandTo(const int size);

	double	Norm(const double l);
	Tp		Projection(Vector &a);
	Tp		DotProduct(Vector &a);

	void	CopyContent(Vector<Tp> &b, const int num);
	virtual Vector  Conv1D( Vector& kern, const bool bcircular);	// circular convolution in 1D
};

////////////////////////////////////////////
// construction & destruction
template<class Tp>
Vector<Tp>::Vector(Vector<Tp>& a){
	mData = a.getData();
	nSize = a.size();
	bRowVector = a.RowVector();
}

////////////////////////////////////////////
// Initialization
template<class Tp> // not tested
void Vector<Tp>::InitValue(Tp value) {
	mData.setTo(cv::Scalar::all(value));
}


////////////////////////////////////////////
// operators

template<class Tp>
Vector<Tp>& Vector<Tp>::operator=(Vector<Tp>& x){
	//beware of self assignments : x = x
	if (this != &x){
		mData = x.getData();
		nSize = x.size();
		bRowVector = x.RowVector();
	}
	return *this;
}

template<class Tp>
Tp& Vector<Tp>::operator ()(const long index){
	/*
	int nRow = index/mData.cols;
	int nCol = index%mData.cols;
	return mData.at<Tp>(nRow, nCol);
	*/
	//cout << "!!!!" << index << endl;
	if (bRowVector){
		return mData.at<Tp>(0, index);
	}
	else{
		return mData.at<Tp>(index, 0);
	}
	//return mData.at<Tp>(index);
}

template<class Tp>
Tp& Vector<Tp>::operator ()(const int row, const int col){
	return mData.at<Tp>(row, col);
}

template<class Tp>
Vector<Tp> Vector<Tp>::operator+ (Vector<Tp>& x){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData + x.getData());
	return temp;
}


template<class Tp>
Vector<Tp> Vector<Tp>::operator+ (const Tp y){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData + y);
	return temp;
}


template<class Tp>
Vector<Tp> Vector<Tp>::operator - (Vector<Tp>& x){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData - x.getData());
	return temp;
}


template<class Tp>
Vector<Tp> Vector<Tp>::operator - (const Tp y){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData - y);
	return temp;
}

// Commented out for testing
/*
template<class Tp>
Vector<Tp> Vector<Tp>::operator * (const Matrix<Tp>& x){
	Vector<Tp> temp;
	temp.setData(mData*x.data);
	return temp;
}
*/

template<class Tp>
Vector<Tp> Vector<Tp>::operator * (const Tp y){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData * y);
	return temp;
}

template<class Tp>
Vector<Tp> Vector<Tp>::operator / (const Tp y){
	Vector<Tp> temp;
	if (bRowVector){
		temp.SetRowVector(size());
	}
	else{
		temp.SetColVector(size());
	}
	temp.setData(mData * (1.0 / (double)y));
	return temp;
}

/*
template<class Tp>
void Vector<Tp>::operator<<(std::ofstream &strm){
	strm << this.mData << endl;
}
*/


////////////////////////////////////////////
// operations
template<class Tp>
Tp Vector<Tp>::Data(const long index){
	if (bRowVector){
		return mData.at<Tp>(0, index);
	}else{
		return mData.at<Tp>(index, 0);
	}
}

template<class Tp>
int Vector<Tp>::rows(){
	return mData.rows;
}

template<class Tp>
int Vector<Tp>::cols(){
	return mData.cols;
}

template<class Tp>
void Vector<Tp>::SetDimension(const int row, const int col){
	if (row == 1){
		SetRowVector(col);
	}else{
		SetColVector(row);
	}
}

template<class Tp>
void Vector<Tp>::Zero(){
	mData = cv::Scalar(0);
}

template<class Tp>
cv::Mat Vector<Tp>::getData(){
	return mData.clone();
}

template<class Tp>
void Vector<Tp>::setData(cv::Mat mNewData){
	mData = mNewData.clone();
}

template<class Tp>
void Vector<Tp>::SetRowVector(const int size){
	nSize = size;
	bRowVector = true;
	mData.create(cv::Size(nSize, 1), cv::DataType<Tp>::type);
}

template<class Tp>
void Vector<Tp>::SetColVector(const int size){
	nSize = size;
	bRowVector = false;
	mData.create(cv::Size(1, nSize), cv::DataType<Tp>::type);
}

template<class Tp>
bool Vector<Tp>::ExpandTo(const int size){
	if (size <= nSize){
		return false;
	}else{
		nSize = size;
	}

	cv::Mat mExpandedData;
	if (bRowVector){
		mExpandedData.create(cv::Size(nSize, 1), cv::DataType<Tp>::type);
	}else{
		mExpandedData.create(cv::Size(1, nSize), cv::DataType<Tp>::type);
	}
	
	mData.copyTo(mExpandedData(cv::Rect(0, 0, mData.cols, mData.rows)));
	mData = mExpandedData;
}

template<class Tp>
bool Vector<Tp>::RowVector(void){
	return bRowVector;
}

template<class Tp>
bool Vector<Tp>::ColVector(void){
	return !bRowVector;
}

template<class Tp>
Tp Vector<Tp>::Projection(Vector<Tp> &a){
	Tp sum = 0;
	for (int i = 0; i<nSize; i++)
		sum = sum + Data(i)*a.Data(i);
	return sum;
}

template<class Tp>
Tp Vector<Tp>::DotProduct(Vector<Tp> &a){
	Tp sum = 0;
	for (int i = 0; i<nSize; i++)
		sum = sum + Data(i)*a.Data(i);
	return sum;
}

template<class Tp>
double Vector<Tp>::Norm(const double l)
{
	if (l == 1){
		return cv::norm(mData, cv::NORM_L1);
	}else if (l == 2){
		return cv::norm(mData, cv::NORM_L2);
	}

	int		i;
	double	sum = 0.0;

	for (i = 0; i<nSize; i++)
		sum = sum + pow(fabs((double)Data(i)), l);
	sum = pow(sum, 1.0 / l);
	return sum;
}

template<class Tp>
void Vector<Tp>::CopyContent(Vector<Tp> &b, const int num)
{
	cv::Mat mNewContent = b.getData()(cv::Rect(0, 0, num, 1));
	if (b.bRowVector){
		mNewContent.copyTo(mData(cv::Rect(0, 0, num, 1)));
	}
	else{
		mNewContent.copyTo(mData(cv::Rect(0, 0, 1, num)));
	}
}

// circular convolution in 1D
template<class Tp>
Vector<Tp> Vector<Tp>::Conv1D(Vector& kern, const bool bcircular)
{
	Vector<Tp>	I;
	long		i, k, m, nsize, nsize_kernel;

	I.SetDimension(rows(), cols());
	I.Zero();

	nsize = size();
	nsize_kernel = kern.size();

	if (bcircular)
	{
		// circular convolution
		for (i = 0; i<nsize; i++)
		{
			for (k = 0; k<nsize_kernel; k++)
			{
				m = (i - (k - (nsize_kernel - 1) / 2) + nsize * 100) % nsize;

				I(i) = I(i) + Data(m)*(kern.Data(k));
			}
		}
	}
	else
	{
		for (i = 0; i<nsize; i++)
		{
			for (k = 0; k<nsize_kernel; k++)
			{
				m = i - (k - (nsize_kernel - 1) / 2);

				if (m >= 0 && m<nsize)
					I(i) = I(i) + Data(m)*(kern.Data(k));
			}
		}
	}

	return I;
}

#endif